## This is the R script for generating the figures and tables
## in both the paper and the appendices


## load libraries --------------------------------------------------------------
## function to load/install multiple libraries
## from: https://gist.github.com/stevenworthington/3178163
ipak <- function(pkg){
  new.pkg <- pkg[!(pkg %in% installed.packages()[, "Package"])]
  if (length(new.pkg)) 
    install.packages(new.pkg, dependencies = TRUE)
  sapply(pkg, require, character.only = TRUE)
}

packages <- c("ggplot2",
"grid",
"xtable",
"gridExtra",
"stargazer",
"plyr",
"foreign",
"countrycode",
"systemfit",
"DataCombine",
"MASS",
"dotwhisker",
"broom",
"dplyr",
"pastecs",
"ggpubr")

ipak(packages)

## multi graph function
grid_arrange_shared_legend <- function(..., ncol = length(list(...)), nrow = 1, position = c("bottom", "right")) {
  
  plots <- list(...)
  position <- match.arg(position)
  g <- ggplotGrob(plots[[1]] + theme(legend.position = position))$grobs
  legend <- g[[which(sapply(g, function(x) x$name) == "guide-box")]]
  lheight <- sum(legend$height)
  lwidth <- sum(legend$width)
  gl <- lapply(plots, function(x) x + theme(legend.position="none"))
  gl <- c(gl, ncol = ncol, nrow = nrow)
  
  combined <- switch(position,
                     "bottom" = arrangeGrob(do.call(arrangeGrob, gl),
                                            legend,
                                            ncol = 1,
                                            heights = unit.c(unit(1, "npc") - lheight, lheight)),
                     "right" = arrangeGrob(do.call(arrangeGrob, gl),
                                           legend,
                                           ncol = 2,
                                           widths = unit.c(unit(1, "npc") - lwidth, lwidth)))
  
  grid.newpage()
  grid.draw(combined)
  
  # return gtable invisibly
  invisible(combined)
  
}


###################################
## Equations for plotting
###################################

## function to plot coefs -----------------------------------------------------
## ----------------------------------------------------------------------------

coef_fun <- function(termname=termname_f_d_l, model=fitsur.f.d.l, labelname=c("Cooperation", "Conflict")){
  coef_coop <- summary(model)$coefficients[grepl("coop_",rownames(summary(model)$coefficients))&!grepl("year",rownames(summary(model)$coefficients)),] 
  coef_conf <- summary(model)$coefficients[grepl("conf_",rownames(summary(model)$coefficients))&!grepl("year",rownames(summary(model)$coefficients)),]
  coef_coop <- coef_coop[-1,] %>% data.frame() %>% tibble::rownames_to_column("term") %>% dplyr::rename(estimate = Estimate, std.error = Std..Error) %>% mutate(model=0, term=termname)
  coef_conf <- coef_conf[-1,] %>% data.frame() %>% tibble::rownames_to_column("term") %>% dplyr::rename(estimate = Estimate, std.error = Std..Error) %>% mutate(model=1, term=termname)
  coef_both <- rbind(coef_coop,coef_conf)
  ind <- grepl("KOF",coef_both$term) | grepl("Major",coef_both$term) 
  coef_both1 <- coef_both[ind,]; coef_both2 <- coef_both[!ind,]
  coefplot1 <- dwplot(coef_both1,
                      vline = geom_vline(xintercept = 0, colour = "grey60", linetype = 2), # plot line at zero _behind_ coefs
                      dot_args = list(aes(shape = model)),
                      whisker_args = list(aes(linetype = model))) +
    theme_bw(base_size = 10) +     xlab("") + ylab("") +
    labs(title="Key Independent Variables") +
    theme(legend.position="bottom",plot.title = element_text(face="bold"), axis.title.y = element_blank(), panel.border = element_blank(), panel.grid.major = element_blank(), panel.grid.minor = element_blank(), axis.line = element_line(colour = "black"))+
    scale_colour_grey(start = 0, end = 0,breaks = c(0, 1),
                      labels = labelname)+
    scale_shape_discrete(breaks = c(0, 1),
                         labels = labelname)
  coefplot2 <- dwplot(coef_both2,
                      vline = geom_vline(xintercept = 0, colour = "grey60", linetype = 2), # plot line at zero _behind_ coefs
                      dot_args = list(aes(shape = model)),
                      whisker_args = list(aes(linetype = model))) +
    theme_bw(base_size = 10) +    xlab("") + ylab("") +
    labs(title="Control Variables") +
    theme(legend.position="bottom",plot.title = element_text(face="bold"), axis.title.y = element_blank(), panel.border = element_blank(), panel.grid.major = element_blank(), panel.grid.minor = element_blank(), axis.line = element_line(colour = "black"))+
    scale_colour_grey(start = 0, end = 0,breaks = c(0, 1),
                      labels = labelname)+
    scale_shape_discrete(breaks = c(0, 1),
                         labels = labelname)
  
  grid_arrange_shared_legend(coefplot1, coefplot2, nrow=1)
}

## function for plotting marginal effects -------------------------------------
## input: data, model
##        para: the major partner type (variable) to compare against the benchmark (both nonmajor)
##        actor: 1 or 2 (selecting conditional on self (1) or target (2)'s KOF)
## ----------------------------------------------------------------------------

margin_fun <- function(data=final.data.dyad, scale01=FALSE, model=fitsur.f.d, 
                       para="major.sym", actor=1, len=100, 
                       labelname=c("Cooperation", "Conflict")){
  if (scale01==TRUE) data = data %>% mutate(KOFEcGIdf.x.0=KOFEcGIdf.x, KOFEcGIdf.y.0=KOFEcGIdf.y)
  glob = seq(quantile(data$KOFEcGIdf.x.0, probs=.25, na.rm=TRUE), quantile(data$KOFEcGIdf.x.0, probs=.75, na.rm=TRUE), length=len)
  set.seed(11)
  coef_sim = as.data.frame(mvrnorm(n=1000, mu=coef(model), Sigma=vcov(model)))[grepl(para, names(coef(model)))&grepl("coop_", names(coef(model)))][c(1,1+actor)]
  coop_plot_dat = data.frame(KOF=glob, margin.mean=NA, margin.se=NA)
  for (i in 1:len){
    coop_plot_dat$margin.mean[i] = mean(exp(coef_sim[,1]+glob[i]*coef_sim[,2])-1)
    coop_plot_dat$margin.se[i] = sd(exp(coef_sim[,1]+glob[i]*coef_sim[,2])-1)
  }
  glob = seq(quantile(data$KOFEcGIdf.x.0, probs=.25, na.rm=TRUE), quantile(data$KOFEcGIdf.x.0, probs=.75, na.rm=TRUE), length=len)
  set.seed(11)
  coef_sim = as.data.frame(mvrnorm(n=1000, mu=coef(model), Sigma=vcov(model)))[grepl(para, names(coef(model)))&grepl("conf", names(coef(model)))][c(1,1+actor)]
  conf_plot_dat = data.frame(KOF=glob, margin.mean=NA, margin.se=NA)
  for (i in 1:len){
    conf_plot_dat$margin.mean[i] = mean(exp(coef_sim[,1]+glob[i]*coef_sim[,2])-1)
    conf_plot_dat$margin.se[i] = sd(exp(coef_sim[,1]+glob[i]*coef_sim[,2])-1)
  }
  plot_dat = rbind(coop_plot_dat, conf_plot_dat) %>% mutate(type=rep(labelname, each=len))
  margin_plot = ggplot(plot_dat, aes(x=KOF, y=100*margin.mean, group=type)) + 
    geom_line(aes(linetype=type),size=1) + geom_ribbon(aes(ymin = 100*(margin.mean-1.96*margin.se), ymax = 100*(margin.mean+1.96*margin.se)), fill = "grey70",alpha = 0.5)+
    theme_bw(base_size = 15) + theme(plot.title = element_text(size = 12),panel.border = element_blank(), panel.grid.major = element_blank(), panel.grid.minor = element_blank(), axis.line = element_line(colour = "black"))+
    xlab("KOF") + ylab(expression("%"~Delta~"Goldstein-weighted events")) + theme(legend.position = "bottom")+
    geom_hline(yintercept = 0, colour = "grey60", linetype = 2)+
    scale_linetype_manual(name="", values=c("dotted","solid"))
  return(margin_plot)
}


## function for plotting marginal effects -------------------------------------
## input: data, model
##        para: the major partner type (variable) to compare against the benchmark (both nonmajor)
##        actor: 1 or 2 (selecting conditional on self (1) or target (2)'s KOF)
## ----------------------------------------------------------------------------

margin_fun2 <- function(data=final.data.dyad, scale01=FALSE, model=fitsur.f.d, 
                       para="KOFEcGIdf.x.0", actor=1, len=100, 
                       labelname=c("Cooperation", "Conflict")){
  if (scale01==TRUE) data = data %>% mutate(KOFEcGIdf.x.0=KOFEcGIdf.x, KOFEcGIdf.y.0=KOFEcGIdf.y)
  set.seed(11)
  glob=apply(data[,para],2, sd, na.rm=TRUE)
  coef_sim = as.data.frame(mvrnorm(n=1000, mu=coef(model), Sigma=vcov(model)))[grepl(para, names(coef(model)))&grepl("coop_", names(coef(model)))]
  margin_sim = glob*data.frame(coef_sim[,1],coef_sim[,1]+coef_sim[,2:4])
  quantilemin = function(x) quantile(x,probs=.025,na.rm=TRUE)
  quantilemax = function(x) quantile(x,probs=.975,na.rm=TRUE)
  coop_plot_dat = data.frame(margin.mean=100*(exp(as.vector(colMeans(margin_sim)))-1),
                             margin.min=100*(exp(apply(margin_sim, 2, quantilemin))-1),
                             margin.max=100*(exp(apply(margin_sim, 2, quantilemax))-1),
                             type=factor(c("NMTP both","MTP both","MTP self","MTP target"),
                                         levels = c("NMTP both","MTP self","MTP target","MTP both")))
  set.seed(11)
  glob=apply(data[,para],2, sd, na.rm=TRUE)
  coef_sim = as.data.frame(mvrnorm(n=1000, mu=coef(model), Sigma=vcov(model)))[grepl(para, names(coef(model)))&grepl("conf_", names(coef(model)))]
  margin_sim = glob*data.frame(coef_sim[,1],coef_sim[,1]+coef_sim[,2:4])
  quantilemin = function(x) quantile(x,probs=.025,na.rm=TRUE)
  quantilemax = function(x) quantile(x,probs=.975,na.rm=TRUE)
  conf_plot_dat = data.frame(margin.mean=100*(exp(as.vector(colMeans(margin_sim)))-1),
                             margin.min=100*(exp(apply(margin_sim, 2, quantilemin))-1),
                             margin.max=100*(exp(apply(margin_sim, 2, quantilemax))-1),
                             type=factor(c("NMTP both","MTP both","MTP self","MTP target"),
                                         levels = c("NMTP both","MTP self","MTP target","MTP both")))
  plot_dat = rbind(coop_plot_dat, conf_plot_dat) %>%
    mutate(eventtype=rep(c("Cooperation", "Conflict"),each=4))
  margin_plot = ggplot(plot_dat, aes(x=type, y=margin.mean, color=eventtype))+ 
    geom_pointrange(aes(ymin = margin.min, ymax = margin.max, shape=eventtype), size=1,position=position_dodge(0.5))+
    theme_bw(base_size = 15) + theme(plot.title = element_text(size = 12),panel.border = element_blank(), panel.grid.major = element_blank(), panel.grid.minor = element_blank(), axis.line = element_line(colour = "black"))+
    xlab("") + ylab(expression("%"~Delta~"Goldstein-weighted events")) + 
    theme(legend.position = "bottom",legend.title=element_blank())+
    geom_hline(yintercept = 0, colour = "red", linetype = 2)+
    scale_linetype_manual(values=c("dashed","solid"))+
    scale_color_manual(values=c('#999999','#E69F00'))
  return(margin_plot)
}


## function for generating prediction plots -----------------------------------
## input: time: select which year to generate predictions
##        traderange: specify the the levels of trade volume
## ----------------------------------------------------------------------------

predict_fun <- function(data=final.data.dyad, scale01=FALSE,time=2005, felevels=1995:2012,
                        model=fitsur.f.d, len=100, se=TRUE, 
                        labelname=c("Cooperation", "Conflict")) {
  if (scale01==TRUE) data = data %>% mutate(KOFEcGIdf.x.0=KOFEcGIdf.x, KOFEcGIdf.y.0=KOFEcGIdf.y)
  temp <- c()
  for (i in 1:length(names(data))){
    temp[i]= sum(grep(names(data)[i],names(coef(model))),na.rm = TRUE)!=0
  }
  new.dat <- apply(data[,names(data)[temp]], 2, median, na.rm=TRUE)
  new.dat <- as.data.frame(t(new.dat)) %>% mutate(coopTotal=mean(data$coopTotal,na.rm=TRUE), 
                                                  hostTotal=mean(data$hostTotal,na.rm=TRUE))
  if (scale01==TRUE) {
    new.dat.ctr <- subset(new.dat,select=-c(major.sym, major.asym1, major.asym2, KOFEcGIdf.x))
  } else{
    new.dat.ctr <- subset(new.dat,select=-c(major.sym, major.asym1, major.asym2, KOFEcGIdf.x.0))
  }
  new.dat.iv = data.frame(KOFEcGIdf.x.0 = seq(quantile(data$KOFEcGIdf.x.0, probs=.25, na.rm=TRUE), quantile(data$KOFEcGIdf.x.0, probs=.75, na.rm=TRUE), length=len)) 
  new.dat.1 = new.dat.iv %>% mutate(major.sym=0, major.asym1=0, major.asym2=0, PartnerType="NMTP both")
  new.dat.2 = new.dat.iv %>% mutate(major.sym=1, major.asym1=0, major.asym2=0, PartnerType="MTP both")
  new.dat.3 = new.dat.iv %>% mutate(major.sym=0, major.asym1=1, major.asym2=0, PartnerType="MTP self")
  new.dat.4 = new.dat.iv %>% mutate(major.sym=0, major.asym1=0, major.asym2=1, PartnerType="MTP target")
  new.dat = bind_rows(new.dat.1,new.dat.2,new.dat.3,new.dat.4) %>% 
    data.frame(., new.dat.ctr) %>% mutate(year=factor(time, levels=felevels)) ## set the year as 2005
  if (scale01==TRUE) new.dat = new.dat %>% mutate(KOFEcGIdf.x=KOFEcGIdf.x.0)
  plot_datf <- cbind(predict(model,newdata=new.dat,se.fit=se),
                     KOF=new.dat$KOFEcGIdf.x.0,
                     PartnerType=new.dat$PartnerType)
  allin1_interact_coop <- ggplot(plot_datf, aes(x=KOF, y=exp(coop.pred)-logdelta, group=PartnerType))+
    geom_line(aes(linetype=PartnerType),size=1) + scale_linetype_manual(name="",values=c("solid", "dotdash", "dashed", "dotted"))+ 
    geom_ribbon(aes(ymin = exp(coop.pred-1.96*coop.se.fit)-logdelta, ymax = exp(coop.pred+1.96*coop.se.fit)-logdelta), fill = "grey70",alpha = 0.5)+
    theme_bw(base_size = 15) + theme(panel.border = element_blank(), legend.key.width = unit(1.5,"cm"), panel.grid.major = element_blank(), panel.grid.minor = element_blank(), axis.line = element_line(colour = "black"))+
    xlab("KOF") + ylab("Goldstein-weighted events") +
    theme(legend.position = "bottom")
  allin1_interact_conf <- ggplot(plot_datf, aes(x=KOF, y=exp(conf.pred)-logdelta, group=PartnerType))+
    geom_line(aes(linetype=PartnerType),size=1) + scale_linetype_manual(name="",values=c("solid", "dotdash", "dashed", "dotted"))+ 
    geom_ribbon(aes(ymin = exp(conf.pred-1.96*conf.se.fit)-logdelta, ymax = exp(conf.pred+1.96*conf.se.fit)-logdelta), fill = "grey70",alpha = 0.5)+
    theme_bw(base_size = 15) + theme(panel.border = element_blank(), legend.key.width = unit(1.5,"cm"),panel.grid.major = element_blank(), panel.grid.minor = element_blank(), axis.line = element_line(colour = "black"))+
    xlab("KOF") + ylab("Goldstein-weighted events") +
    theme(legend.position = "bottom")
  allin1_interact <- grid_arrange_shared_legend(allin1_interact_coop+labs(title=paste("(a)", labelname[1])), 
                                                allin1_interact_conf+labs(title=paste("(b)", labelname[2])), nrow=1)
  return(allin1_interact)
}

## function for generating Latex tables
## modelname: specify as the directory for the folder in the appendix
## formfake: specify one of the formulas
## termname: output for the table's covariates name
## lablename: output for the table's column name
latex_fun <- function(model=fitsur.f.d.l, modelname="DYLDVModel", data=final.data.dyad, 
                      formfake=coop.f.d.l, termname=termname_f_d_l, fontsize="normalsize",
                      labelname=c("Cooperation", "Conflict"), surols="OLS", 
                      reportfe=FALSE, fename="Year",
                      modelcap="lagged events model"){
  vec.coef <- list()
  vec.se <- list()
  vec.t <- list()
  vec.ci <- list()
  for (i in 1:2){
    table <- summary(model$eq[[i]])$coefficients
    vec.coef[[i]] <- as.matrix(table[,1]) 
    vec.se[[i]] <- as.matrix(table[,2])
    vec.ci[[i]] <- as.matrix(cbind(vec.coef[[i]]-1.96*vec.se[[i]],vec.coef[[i]]+1.96*vec.se[[i]])) 
    vec.t[[i]] <- vec.coef[[i]]/vec.se[[i]]
  }
  ## now fake lm objects
  p1 <- lm(formfake, data)
  if (reportfe==TRUE){
    term.cut <- grep("year",rownames(summary(model)$coefficients))[1]-2 # index for the last term before year
    vec.names <- c(termname[1:term.cut],paste0(fename,":",1996:2012),termname[-c(1:term.cut)], "Constant")
    omityear <- NULL
    outfile <-paste0("AppendixMaterials/", modelname,"/", modelname, "wFE.tex")
  } else {
    vec.names <- c(termname, "Constant")
    omityear <- "year"
    outfile <-paste0("AppendixMaterials/", modelname,"/", modelname, ".tex")
  }
  stargazer(p1,p1, title=paste0(surols, " coefficients and 95 percent confidence bounds for ", modelcap, ", 1995-2012"), 
            single.row = TRUE, omit.stat="all", dep.var.caption="", dep.var.labels  = "Dependent Variables",
            column.labels=labelname, #dep.var.labels= labelname, 
            covariate.labels=vec.names, omit=omityear,
            font.size = fontsize,
            no.space=TRUE,
            model.numbers = FALSE,
            digits=2,
            coef = vec.coef,
            ci.custom = vec.ci,
            t = vec.t, 
            add.lines = list(c("Observations", sum(summary(model$eq[[1]])$df),sum(summary(model$eq[[2]])$df)),
                             c("R$^{2}$",round(c(sum(summary(model$eq[[1]])$r.squared),sum(summary(model$eq[[2]])$r.squared)),3)),
                             c("Adjusted R$^{2}$",round(c(sum(summary(model$eq[[1]])$adj.r.squared),sum(summary(model$eq[[2]])$adj.r.squared)),3))),
            out.header = FALSE, out=outfile
  )
}


################################
## Main Model: Dyadic LDV
################################
dir.create("./AppendixMaterials")
dir.create("./AppendixMaterials/DYLDVModel")

## coef plot -----------------------------------------------------------------
## term names for DY LDV model
termname_f_d_l <- c("Coop Level","Conf Level","KOF 1","KOF 2",
                    "Mutual Major Trade Partners", "1 Is Major Partner of 2", "2 Is Major Partner of 1", 
                    "Trade/GDP 1", "Trade/GDP 2", "Trade Balance",
                    "CINC Ratio", "log Distance", "log Population 1",
                    "log Population 2", "log GDP 1", "log GDP 2",
                    "Democracy 1", "Democracy 2",
                    "Alliance","Contiguity",
                    "KOF 1:Mutual Major Trade Partners",
                    "KOF 1:1 Is Major Partner of 2", "KOF 1:2 Is Major Partner of 1",
                    "KOF 2:Mutual Major Trade Partners",
                    "KOF 2:1 Is Major Partner of 2", "KOF 2:2 Is Major Partner of 1",
                    "Trade/GDP 1:Trade/GDP 2", "Democracy 1:Democracy 2")

coefplot_f_d_l <- coef_fun(termname=termname_f_d_l, model=fitsur.f.d.l)

ggsave(file = "AppendixMaterials/DYLDVModel/coefplot_f_d_l.pdf", coefplot_f_d_l, width=10, height=6)
ggsave(file = "Figures/coefplot_f_d_l.pdf", coefplot_f_d_l, width=10, height=6)


## Marginal Impact Plots ------------------------------------------------------

margin_plot_bothmajor_f_d_l <- margin_fun(model=fitsur.f.d.l,para="major.sym")
margin_plot_selfmajor_f_d_l <- margin_fun(model=fitsur.f.d.l,para="major.asym1")
margin_plot_targetmajor_f_d_l <- margin_fun(model=fitsur.f.d.l,para="major.asym2")
margin_plot_f_d_l <- grid_arrange_shared_legend(margin_plot_bothmajor_f_d_l+labs(title="(a) MTP both"), 
                                                margin_plot_selfmajor_f_d_l +labs(title="(b) MTP self"), 
                                                margin_plot_targetmajor_f_d_l+labs(title="(c) MTP target"), nrow=1)
ggsave(file = "AppendixMaterials/DYLDVModel/marginplot_f_d_l.pdf", margin_plot_f_d_l, width=12, height=6)
ggsave(file = "Figures/marginplot_f_d_l.pdf", margin_plot_f_d_l, width=12, height=6)

margin_plot_kof_f_d_l <- margin_fun2(model=fitsur.f.d.l)
ggsave(file = "AppendixMaterials/DYLDVModel/marginplot_kof_f_d_l.pdf", margin_plot_kof_f_d_l, width=12, height=6)
ggsave(file = "Figures/marginplot_kof_f_d_l.pdf", margin_plot_kof_f_d_l, width=12, height=6)

margin_plot_allin1_f_d_l <- grid.arrange(margin_plot_f_d_l,
                                         margin_plot_kof_f_d_l+labs(title="(d) KOF"),
                                         nrow=2)
ggsave(file = "AppendixMaterials/DYLDVModel/marginplot_allin1_f_d_l.pdf", margin_plot_allin1_f_d_l, width=12, height=10)
ggsave(file = "Figures/marginplot_allin1_f_d_l.pdf", margin_plot_allin1_f_d_l, width=12, height=10)



## Prediction Plot ------------------------------------------------------------

predict_plot_f_d_l <- predict_fun(model=fitsur.f.d.l)
ggsave(file = "AppendixMaterials/DYLDVModel/predictplot_f_d_l.pdf", predict_plot_f_d_l, width=10, height=6)
ggsave(file = "Figures/predictplot_f_d_l.pdf", predict_plot_f_d_l, width=10, height=6)


## Main Table: DY LDV ---------------------------------------------------------
latex_fun(model=fitsur.f.d.l, modelname="DYLDVModel", modelcap="lagged events model", data=final.data.dyad, formfake=coop.f.d.l, termname=termname_f_d_l, labelname=c("Cooperation", "Conflict"))


## ------------------------------------------------------------------
## SUR MODEL
dir.create("./AppendixMaterials/DYLDVModelsur")

## coef plot -----------------------------------------------------------------
## term names for DY LDV model
termname_sur_f_d_l <- c("Coop Level","Conf Level","KOF 1","KOF 2",
                    "Mutual Major Trade Partners", "1 Is Major Partner of 2", "2 Is Major Partner of 1", 
                    "Trade/GDP 1", "Trade/GDP 2", "Trade Balance",
                    "CINC Ratio", "log Distance", "log Population 1",
                    "log Population 2", "log GDP 1", "log GDP 2",
                    "Democracy 1", "Democracy 2",
                    "Alliance","Contiguity",
                    "KOF 1:Mutual Major Trade Partners",
                    "KOF 1:1 Is Major Partner of 2", "KOF 1:2 Is Major Partner of 1",
                    "KOF 2:Mutual Major Trade Partners",
                    "KOF 2:1 Is Major Partner of 2", "KOF 2:2 Is Major Partner of 1",
                    "Trade/GDP 1:Trade/GDP 2", "Democracy 1:Democracy 2")

coefplot_sur_f_d_l <- coef_fun(termname=termname_sur_f_d_l, model=sur.f.d.l)

ggsave(file = "AppendixMaterials/DYLDVModelsur/coefplot_sur_f_d_l.pdf", coefplot_sur_f_d_l, width=10, height=6)


## Marginal Impact Plots ------------------------------------------------------

margin_plot_bothmajor_sur_f_d_l <- margin_fun(model=sur.f.d.l,para="major.sym")
margin_plot_selfmajor_sur_f_d_l <- margin_fun(model=sur.f.d.l,para="major.asym1")
margin_plot_targetmajor_sur_f_d_l <- margin_fun(model=sur.f.d.l,para="major.asym2")
margin_plot_sur_f_d_l <- grid_arrange_shared_legend(margin_plot_bothmajor_sur_f_d_l+labs(title="(a) MTP both"), 
                                                margin_plot_selfmajor_sur_f_d_l +labs(title="(b) MTP self"), 
                                                margin_plot_targetmajor_sur_f_d_l+labs(title="(c) MTP target"), nrow=1)
ggsave(file = "AppendixMaterials/DYLDVModelsur/marginplot_sur_f_d_l.pdf", margin_plot_sur_f_d_l, width=12, height=6)

margin_plot_kof_sur_f_d_l <- margin_fun2(model=sur.f.d.l)
ggsave(file = "AppendixMaterials/DYLDVModelsur/marginplot_kof_sur_f_d_l.pdf", margin_plot_kof_sur_f_d_l, width=12, height=6)

margin_plot_allin1_sur_f_d_l <- grid.arrange(margin_plot_sur_f_d_l,
                                         margin_plot_kof_sur_f_d_l+labs(title="(d) KOF"),
                                         nrow=2)
ggsave(file = "AppendixMaterials/DYLDVModelsur/marginplot_allin1_sur_f_d_l.pdf", margin_plot_allin1_sur_f_d_l, width=12, height=10)


## Prediction Plot ------------------------------------------------------------

predict_plot_sur_f_d_l <- predict_fun(model=sur.f.d.l)
ggsave(file = "AppendixMaterials/DYLDVModelsur/predictplot_sur_f_d_l.pdf", predict_plot_sur_f_d_l, width=10, height=6)
## table
latex_fun(model=sur.f.d.l, surols = "SUR", modelname="DYLDVModelsur", modelcap="lagged events model", data=final.data.dyad, formfake=coop.f.d.l, termname=termname_sur_f_d_l, labelname=c("Cooperation", "Conflict"))



#################################################
################################################
## Models in the Appendix
##################################################
################################################


###########################
## 1. DY Year FE
###########################
dir.create("./AppendixMaterials/DYFEModel")

## Coef Plot
termname_f_d <- c("KOF 1","KOF 2","Mutual Major Trade Partners",
                  "1 Is Major Partner of 2", "2 Is Major Partner of 1", "Trade/GDP 1", "Trade/GDP 2",
                  "Trade Balance",
                  "CINC Ratio", "log Distance", "log Population 1",
                  "log Population 2", "log GDP PC 1", "log GDP PC 2",
                  "Democracy 1", "Democracy 2", "Alliance","Contiguity", "KOF 1:Mutual Major Trade Partners",
                  "KOF 1:1 Is Major Partner of 2", "KOF 1:2 Is Major Partner of 1",
                  "KOF 2:Mutual Major Trade Partners",
                  "KOF 2:1 Is Major Partner of 2", "KOF 2:2 Is Major Partner of 1",
                  "Trade/GDP 1:Trade/GDP 2", "Democracy 1:Democracy 2")

coefplot_f_d <- coef_fun(termname=termname_f_d, model=fitsur.f.d)

ggsave(file = "AppendixMaterials/DYFEModel/coefplot_f_d.pdf", coefplot_f_d, width=10, height=6)


## Marginal Effects plot
margin_plot_bothmajor_f_d <- margin_fun(model=fitsur.f.d,para="major.sym")
margin_plot_selfmajor_f_d <- margin_fun(model=fitsur.f.d,para="major.asym1")
margin_plot_targetmajor_f_d <- margin_fun(model=fitsur.f.d,para="major.asym2")
margin_plot_f_d <- grid_arrange_shared_legend(margin_plot_bothmajor_f_d+labs(title="(a) MTP both"), 
                                              margin_plot_selfmajor_f_d+labs(title="(b) MTP self"), 
                                              margin_plot_targetmajor_f_d+labs(title="(c) MTP target"), nrow=1)
ggsave(file = "AppendixMaterials/DYFEModel/marginplot_f_d.pdf", margin_plot_f_d, width=12, height=6)

margin_plot_kof_f_d <- margin_fun2(model=fitsur.f.d)
ggsave(file = "AppendixMaterials/DYFEModel/marginplot_kof_f_d.pdf", margin_plot_kof_f_d, width=12, height=6)

margin_plot_allin1_f_d <- grid.arrange(margin_plot_f_d,
                                         margin_plot_kof_f_d+labs(title="(d) KOF"),
                                         nrow=2)
ggsave(file = "AppendixMaterials/DYFEModel/marginplot_allin1_f_d.pdf", margin_plot_allin1_f_d, width=12, height=10)


## Predict Plot
predict_plot_f_d <- predict_fun(model=fitsur.f.d)
ggsave(file = "AppendixMaterials/DYFEModel/predictplot_f_d.pdf", predict_plot_f_d, width=10, height=6)


## Table
latex_fun(model=fitsur.f.d, modelname="DYFEModel", modelcap="model with year fixed effects", data=final.data.dyad, formfake=coop.f.d, termname=termname_f_d, labelname=c("Cooperation", "Conflict"))

## Table with fixed effects
latex_fun(model=fitsur.f.d, reportfe=TRUE, fontsize = "footnotesize", modelname="DYFEModel", modelcap="model with year fixed effects", data=final.data.dyad, formfake=coop.f.d, termname=termname_f_d, labelname=c("Cooperation", "Conflict"))


###########################
## 2. DY LDV three way interaction
###########################
dir.create("./AppendixMaterials/DYLDVIModel")

## Coef Plot
termname_f_d_l_i <- c("Coop Level", "Conf Level",  
                      "KOF 1","KOF 2","Mutual Major Trade Partners",
                      "1 Is Major Partner of 2", "2 Is Major Partner of 1", "Trade/GDP 1", "Trade/GDP 2",
                      "Trade Balance", "CINC Ratio", "log Distance", 
                      "log Population 1","log Population 2", "log GDP PC 1", "log GDP PC 2",
                      "Democracy 1", "Democracy 2", "Alliance","Contiguity", 
                      "KOF 1:KOF2", "KOF 1:Mutual Major Trade Partners",
                      "KOF 1:1 Is Major Partner of 2", "KOF 1:2 Is Major Partner of 1",
                      "KOF 2:Mutual Major Trade Partners",
                      "KOF 2:1 Is Major Partner of 2", "KOF 2:2 Is Major Partner of 1",
                      "Trade/GDP 1:Trade/GDP 2", "Democracy 1:Democracy 2",
                      "KOF 1:KOF2:Mutual Major Trade Partners",
                      "KOF 1:KOF2:1 Is Major Partner of 2",
                      "KOF 1:KOF2:2 Is Major Partner of 1")

coefplot_f_d_l_i <- coef_fun(termname=termname_f_d_l_i, model=fitsur.f.d.l.i)

ggsave(file = "AppendixMaterials/DYLDVIModel/coefplot_f_d_l_i.pdf", coefplot_f_d_l_i, width=10, height=6)


## Marginal Effects plot
margin_plot_bothmajor_f_d_l_i <- margin_fun(model=fitsur.f.d.l.i,para="major.sym")
margin_plot_selfmajor_f_d_l_i <- margin_fun(model=fitsur.f.d.l.i,para="major.asym1")
margin_plot_targetmajor_f_d_l_i <- margin_fun(model=fitsur.f.d.l.i,para="major.asym2")
margin_plot_f_d_l_i <- grid_arrange_shared_legend(margin_plot_bothmajor_f_d_l_i+labs(title="(a) MTP both"), 
                                                  margin_plot_selfmajor_f_d_l_i+labs(title="(b) MTP self"), 
                                                  margin_plot_targetmajor_f_d_l_i+labs(title="(c) MTP target"), nrow=1)
ggsave(file = "AppendixMaterials/DYLDVIModel/marginplot_f_d_l_i.pdf", margin_plot_f_d_l_i, width=12, height=6)


margin_plot_kof_f_d_l_i <- margin_fun2(model=fitsur.f.d.l.i)
ggsave(file = "AppendixMaterials/DYLDVIModel/marginplot_kof_f_d_l_i.pdf", margin_plot_kof_f_d_l_i, width=12, height=6)

margin_plot_allin1_f_d_l_i <- grid.arrange(margin_plot_f_d_l_i,
                                       margin_plot_kof_f_d_l_i+labs(title="(d) KOF"),
                                       nrow=2)
ggsave(file = "AppendixMaterials/DYLDVIModel/marginplot_allin1_f_d_l_i.pdf", margin_plot_allin1_f_d_l_i, width=12, height=10)


## Predict Plot
predict_plot_f_d_l_i <- predict_fun(model=fitsur.f.d.l.i)
ggsave(file = "AppendixMaterials/DYLDVIModel/predictplot_f_d_l_i.pdf", predict_plot_f_d_l_i, width=10, height=6)

## Table
latex_fun(model=fitsur.f.d.l.i, modelname="DYLDVIModel", modelcap="lagged events model with three-way interaction",data=final.data.dyad, formfake=coop.f.d.l.i, termname=termname_f_d_l_i, labelname=c("Cooperation", "Conflict"))


###########################
## 3. DY LDV relative events
###########################
dir.create("./AppendixMaterials/DYLDVTModel")

## Coef Plot
termname_f_d_l_t <- c("Coop Level","Conf Level","KOF 1","KOF 2",
                      "Mutual Major Trade Partners", "1 Is Major Partner of 2", "2 Is Major Partner of 1", 
                      "Trade/GDP 1", "Trade/GDP 2", "Trade Balance",
                      "CINC Ratio", "log Distance", "log Population 1",
                      "log Population 2", "log GDP PC 1", "log GDP PC 2",
                      "Democracy 1", "Democracy 2", "Alliance","Contiguity", "KOF 1:Mutual Major Trade Partners",
                      "KOF 1:1 Is Major Partner of 2", "KOF 1:2 Is Major Partner of 1",
                      "KOF 2:Mutual Major Trade Partners",
                      "KOF 2:1 Is Major Partner of 2", "KOF 2:2 Is Major Partner of 1",
                      "Trade/GDP 1:Trade/GDP 2", "Democracy 1:Democracy 2")

coefplot_f_d_l_t <- coef_fun(termname=termname_f_d_l_t, model=fitsur.f.d.l.t, labelname = c("Relative","Total"))

ggsave(file = "AppendixMaterials/DYLDVTModel/coefplot_f_d_l_t.pdf", coefplot_f_d_l_t, width=10, height=6)


## Marginal Effects plot
margin_plot_bothmajor_f_d_l_t <- margin_fun(model=fitsur.f.d.l.t,para="major.sym", labelname = c("Relative","Total"))
margin_plot_selfmajor_f_d_l_t <- margin_fun(model=fitsur.f.d.l.t,para="major.asym1", labelname = c("Relative","Total"))
margin_plot_targetmajor_f_d_l_t <- margin_fun(model=fitsur.f.d.l.t,para="major.asym2", labelname = c("Relative","Total"))
margin_plot_f_d_l_t <- grid_arrange_shared_legend(margin_plot_bothmajor_f_d_l_t+labs(title="(a) MTP both"), 
                                                  margin_plot_selfmajor_f_d_l_t+labs(title="(b) MTP self"), 
                                                  margin_plot_targetmajor_f_d_l_t+labs(title="(c) MTP target"), nrow=1)
ggsave(file = "AppendixMaterials/DYLDVTModel/marginplot_f_d_l_t.pdf", margin_plot_f_d_l_t, width=12, height=6)


margin_plot_kof_f_d_l_t <- margin_fun2(model=fitsur.f.d.l.t)
ggsave(file = "AppendixMaterials/DYLDVTModel/marginplot_kof_f_d_l_t.pdf", margin_plot_kof_f_d_l_t, width=12, height=6)

margin_plot_allin1_f_d_l_t <- grid.arrange(margin_plot_f_d_l_t,
                                           margin_plot_kof_f_d_l_t+labs(title="(d) KOF"),
                                           nrow=2)
ggsave(file = "AppendixMaterials/DYLDVTModel/marginplot_allin1_f_d_l_t.pdf", margin_plot_allin1_f_d_l_t, width=12, height=10)



## Predict Plot
predict_plot_f_d_l_t <- predict_fun(model=fitsur.f.d.l.t, labelname = c("Relative","Total"))
ggsave(file = "AppendixMaterials/DYLDVTModel/predictplot_f_d_l_t.pdf", predict_plot_f_d_l_t, width=10, height=6)


## Table
latex_fun(model=fitsur.f.d.l.t, modelname="DYLDVTModel", modelcap="relative/total model", data=final.data.dyad, formfake=rel.f.d.l.t, termname=termname_f_d_l_t, labelname=c("Relative","Total"))



###########################
## 4.DY Year FE without Trade DEP variables
###########################
dir.create("./AppendixMaterials/DYFEWODModel")

## Coef Plot
termname_f_d_wod <- c("KOF 1","KOF 2","Mutual Major Trade Partners",
                      "1 Is Major Partner of 2", "2 Is Major Partner of 1", 
                      "Trade Balance",
                      "CINC Ratio", "log Distance", "log Population 1",
                      "log Population 2", "log GDP PC 1", "log GDP PC 2",
                      "Democracy 1", "Democracy 2", "Alliance","Contiguity", "KOF 1:Mutual Major Trade Partners",
                      "KOF 1:1 Is Major Partner of 2", "KOF 1:2 Is Major Partner of 1",
                      "KOF 2:Mutual Major Trade Partners",
                      "KOF 2:1 Is Major Partner of 2", "KOF 2:2 Is Major Partner of 1",
                      "Democracy 1:Democracy 2")

coefplot_f_d_wod <- coef_fun(termname=termname_f_d_wod, model=fitsur.f.d.wod)

ggsave(file = "AppendixMaterials/DYFEWODModel/coefplot_f_d_wod.pdf", coefplot_f_d_wod, width=10, height=6)


## Marginal Effects plot
margin_plot_bothmajor_f_d_wod <- margin_fun(model=fitsur.f.d.wod,para="major.sym")
margin_plot_selfmajor_f_d_wod <- margin_fun(model=fitsur.f.d.wod,para="major.asym1")
margin_plot_targetmajor_f_d_wod <- margin_fun(model=fitsur.f.d.wod,para="major.asym2")
margin_plot_f_d_wod <- grid_arrange_shared_legend(margin_plot_bothmajor_f_d_wod+labs(title="(a) MTP both"), 
                                                  margin_plot_selfmajor_f_d_wod+labs(title="(b) MTP self"), 
                                                  margin_plot_targetmajor_f_d_wod+labs(title="(c) MTP target"), nrow=1)
ggsave(file = "AppendixMaterials/DYFEWODModel/marginplot_f_d_wod.pdf", margin_plot_f_d_wod, width=12, height=6)

margin_plot_kof_f_d_wod <- margin_fun2(model=fitsur.f.d.wod)
ggsave(file = "AppendixMaterials/DYFEWODModel/marginplot_kof_f_d_wod.pdf", margin_plot_kof_f_d_wod, width=12, height=6)

margin_plot_allin1_f_d_wod <- grid.arrange(margin_plot_f_d_wod,
                                           margin_plot_kof_f_d_wod+labs(title="(d) KOF"),
                                           nrow=2)
ggsave(file = "AppendixMaterials/DYFEWODModel/marginplot_allin1_f_d_wod.pdf", margin_plot_allin1_f_d_wod, width=12, height=10)





## Predict Plot
predict_plot_f_d_wod <- predict_fun(model=fitsur.f.d.wod)
ggsave(file = "AppendixMaterials/DYFEWODModel/predictplot_f_d_wod.pdf", predict_plot_f_d_wod, width=10, height=6)

## Table
latex_fun(model=fitsur.f.d.wod, modelname="DYFEWODModel", modelcap="year FE, omitted dependence model", data=final.data.dyad, formfake=coop.f.d.wod, termname=termname_f_d_wod, labelname=c("Cooperation", "Conflict"))

## Table with fixed effects
latex_fun(model=fitsur.f.d.wod, reportfe=TRUE, fontsize = "footnotesize", modelname="DYFEWODModel", modelcap="year FE, omitted dependence model", data=final.data.dyad, formfake=coop.f.d.wod, termname=termname_f_d_wod, labelname=c("Cooperation", "Conflict"))


################################
## 5. DY LDV without Trade DEP variables
################################
dir.create("./AppendixMaterials/DYLDVWODModel")

## coef plot -----------------------------------------------------------------
## term names for DY LDV model
termname_f_d_l_wod <- c("Coop Level","Conf Level","KOF 1","KOF 2",
                       "Mutual Major Trade Partners", "1 Is Major Partner of 2", "2 Is Major Partner of 1", 
                       "Trade Balance",
                       "CINC Ratio", "log Distance", "log Population 1",
                       "log Population 2", "log GDP PC 1", "log GDP PC 2",
                       "Democracy 1", "Democracy 2", "Alliance","Contiguity", "KOF 1:Mutual Major Trade Partners",
                       "KOF 1:1 Is Major Partner of 2", "KOF 1:2 Is Major Partner of 1",
                       "KOF 2:Mutual Major Trade Partners",
                       "KOF 2:1 Is Major Partner of 2", "KOF 2:2 Is Major Partner of 1",
                       "Democracy 1:Democracy 2")

coefplot_f_d_l_wod <- coef_fun(termname=termname_f_d_l_wod, model=fitsur.f.d.l.wod)

ggsave(file = "AppendixMaterials/DYLDVWODModel/coefplot_f_d_l_wod.pdf", coefplot_f_d_l_wod, width=10, height=6)


## Marginal Impact Plots ------------------------------------------------------

margin_plot_bothmajor_f_d_l_wod <- margin_fun(model=fitsur.f.d.l.wod,para="major.sym")
margin_plot_selfmajor_f_d_l_wod <- margin_fun(model=fitsur.f.d.l.wod,para="major.asym1")
margin_plot_targetmajor_f_d_l_wod <- margin_fun(model=fitsur.f.d.l.wod,para="major.asym2")
margin_plot_f_d_l_wod <- grid_arrange_shared_legend(margin_plot_bothmajor_f_d_l_wod+labs(title="(a) MTP both"), 
                                                   margin_plot_selfmajor_f_d_l_wod +labs(title="(b) MTP self"), 
                                                   margin_plot_targetmajor_f_d_l_wod+labs(title="(c) MTP target"), nrow=1)
ggsave(file = "AppendixMaterials/DYLDVWODModel/marginplot_f_d_l_wod.pdf", margin_plot_f_d_l_wod, width=12, height=6)

margin_plot_kof_f_d_l_wod <- margin_fun2(model=fitsur.f.d.l.wod)
ggsave(file = "AppendixMaterials/DYLDVWODModel/marginplot_kof_f_d_l_wod.pdf", margin_plot_kof_f_d_l_wod, width=12, height=6)

margin_plot_allin1_f_d_l_wod <- grid.arrange(margin_plot_f_d_l_wod,
                                           margin_plot_kof_f_d_l_wod+labs(title="(d) KOF"),
                                           nrow=2)
ggsave(file = "AppendixMaterials/DYLDVWODModel/marginplot_allin1_f_d_l_wod.pdf", margin_plot_allin1_f_d_l_wod, width=12, height=10)



## Prediction Plot ------------------------------------------------------------

predict_plot_f_d_l_wod <- predict_fun(model=fitsur.f.d.l.wod)
ggsave(file = "AppendixMaterials/DYLDVWODModel/predictplot_f_d_l_wod.pdf", predict_plot_f_d_l_wod, width=10, height=6)

## Table
latex_fun(model=fitsur.f.d.l.wod, modelname="DYLDVWODModel", modelcap="lagged events, omitted dependence model", data=final.data.dyad, formfake=coop.f.d.l.wod, termname=termname_f_d_l_wod, labelname=c("Cooperation", "Conflict"))


################################
## 6. DY LDV with KOF scale between 0 and 1
################################
dir.create("./AppendixMaterials/DYLDV01Model")

## coef plot -----------------------------------------------------------------
## term names for DY LDV model
termname_f_d_l_01 <- termname_f_d_l

coefplot_f_d_l_01 <- coef_fun(termname=termname_f_d_l_01, model=fitsur.f.d.l.01)

ggsave(file = "AppendixMaterials/DYLDV01Model/coefplot_f_d_l_01.pdf", coefplot_f_d_l_01, width=10, height=6)


## Marginal Impact Plots ------------------------------------------------------

margin_plot_bothmajor_f_d_l_01 <- margin_fun(model=fitsur.f.d.l.01,para="major.sym", scale01 = TRUE)
margin_plot_selfmajor_f_d_l_01 <- margin_fun(model=fitsur.f.d.l.01,para="major.asym1", scale01 = TRUE)
margin_plot_targetmajor_f_d_l_01 <- margin_fun(model=fitsur.f.d.l.01,para="major.asym2", scale01 = TRUE)
margin_plot_f_d_l_01 <- grid_arrange_shared_legend(margin_plot_bothmajor_f_d_l_01+labs(title="(a) MTP both"), 
                                                   margin_plot_selfmajor_f_d_l_01 +labs(title="(b) MTP self"), 
                                                   margin_plot_targetmajor_f_d_l_01+labs(title="(c) MTP target"), nrow=1)
ggsave(file = "AppendixMaterials/DYLDV01Model/marginplot_f_d_l_01.pdf", margin_plot_f_d_l_01, width=12, height=6)

margin_plot_kof_f_d_l_01 <- margin_fun2(model=fitsur.f.d.l.01, scale01=TRUE)
ggsave(file = "AppendixMaterials/DYLDV01Model/marginplot_kof_f_d_l_01.pdf", margin_plot_kof_f_d_l_01, width=12, height=6)

margin_plot_allin1_f_d_l_01 <- grid.arrange(margin_plot_f_d_l_01,
                                             margin_plot_kof_f_d_l_01+labs(title="(d) KOF"),
                                             nrow=2)
ggsave(file = "AppendixMaterials/DYLDV01Model/marginplot_allin1_f_d_l_01.pdf", margin_plot_allin1_f_d_l_01, width=12, height=10)



## Prediction Plot ------------------------------------------------------------

predict_plot_f_d_l_01 <- predict_fun(model=fitsur.f.d.l.01, scale01=TRUE)
ggsave(file = "AppendixMaterials/DYLDV01Model/predictplot_f_d_l_01.pdf", predict_plot_f_d_l_01, width=10, height=6)

## Table
latex_fun(model=fitsur.f.d.l.01, modelname="DYLDV01Model", data=final.data.dyad, formfake=coop.f.d.l.01, termname=termname_f_d_l_01, labelname=c("Cooperation", "Conflict"))



################################
## 7. DY LDV with pta
################################
dir.create("./AppendixMaterials/DYLDVModelwpta")

## coef plot -----------------------------------------------------------------
## term names for DY LDV model
termname_f_d_l_pta <- c("Coop Level","Conf Level","KOF 1","KOF 2",
                        "Mutual Major Trade Partners", "1 Is Major Partner of 2", "2 Is Major Partner of 1", 
                        "Trade/GDP 1", "Trade/GDP 2", "Trade Balance",
                        "CINC Ratio", "log Distance", "log Population 1",
                        "log Population 2", "log GDP 1", "log GDP 2",
                        "Democracy 1", "Democracy 2",
                        "Alliance","Contiguity",
                        "Dyad PTAs", "log PTAs 1", "log PTAs 2",
                        "KOF 1:Mutual Major Trade Partners",
                        "KOF 1:1 Is Major Partner of 2", "KOF 1:2 Is Major Partner of 1",
                        "KOF 2:Mutual Major Trade Partners",
                        "KOF 2:1 Is Major Partner of 2", "KOF 2:2 Is Major Partner of 1",
                        "Trade/GDP 1:Trade/GDP 2", "Democracy 1:Democracy 2")

coefplot_f_d_l_pta <- coef_fun(termname=termname_f_d_l_pta, model=fitsur.f.d.l.pta)

ggsave(file = "AppendixMaterials/DYLDVModelwpta/coefplot_f_d_l_pta.pdf", coefplot_f_d_l_pta, width=10, height=6)


## Marginal Impact Plots ------------------------------------------------------

margin_plot_bothmajor_f_d_l_pta <- margin_fun(model=fitsur.f.d.l.pta,para="major.sym")
margin_plot_selfmajor_f_d_l_pta <- margin_fun(model=fitsur.f.d.l.pta,para="major.asym1")
margin_plot_targetmajor_f_d_l_pta <- margin_fun(model=fitsur.f.d.l.pta,para="major.asym2")
margin_plot_f_d_l_pta <- grid_arrange_shared_legend(margin_plot_bothmajor_f_d_l_pta+labs(title="(a) MTP both"), 
                                                   margin_plot_selfmajor_f_d_l_pta +labs(title="(b) MTP self"), 
                                                   margin_plot_targetmajor_f_d_l_pta+labs(title="(c) MTP target"), nrow=1)
ggsave(file = "AppendixMaterials/DYLDVModelwpta/marginplot_f_d_l_pta.pdf", margin_plot_f_d_l_pta, width=12, height=6)

margin_plot_kof_f_d_l_pta <- margin_fun2(model=fitsur.f.d.l.pta)
ggsave(file = "AppendixMaterials/DYLDVModelwpta/marginplot_kof_f_d_l_pta.pdf", margin_plot_kof_f_d_l_pta, width=12, height=6)

margin_plot_allin1_f_d_l_pta <- grid.arrange(margin_plot_f_d_l_pta,
                                            margin_plot_kof_f_d_l_pta+labs(title="(d) KOF"),
                                            nrow=2)
ggsave(file = "AppendixMaterials/DYLDVModelwpta/marginplot_allin1_f_d_l_pta.pdf", margin_plot_allin1_f_d_l_pta, width=12, height=10)



## Prediction Plot ------------------------------------------------------------

predict_plot_f_d_l_pta <- predict_fun(model=fitsur.f.d.l.pta)
ggsave(file = "AppendixMaterials/DYLDVModelwpta/predictplot_f_d_l_pta.pdf", predict_plot_f_d_l_pta, width=10, height=6)

## Table
latex_fun(model=fitsur.f.d.l.pta, modelname="DYLDVModelwpta", data=final.data.dyad, formfake=coop.f.d.l.pta, termname=termname_f_d_l_pta, labelname=c("Cooperation", "Conflict"))


################################
## 8. DY LDV KOF replaced with major partners count
################################
dir.create("./AppendixMaterials/DYLDVModelcount")

## coef plot -----------------------------------------------------------------
## term names for DY LDV model
termname_f_d_l_count <- c("Coop Level","Conf Level","Major Partners Count 1","Major Partners Count 2",
                        "Mutual Major Trade Partners", "1 Is Major Partner of 2", "2 Is Major Partner of 1", 
                        "Trade/GDP 1", "Trade/GDP 2", "Trade Balance",
                        "CINC Ratio", "log Distance", "log Population 1",
                        "log Population 2", "log GDP 1", "log GDP 2",
                        "Democracy 1", "Democracy 2",
                        "Alliance","Contiguity",
                        "Major Partners Count 1:Mutual Major Trade Partners",
                        "Major Partners Count 1:1 Is Major Partner of 2", "Major Partners Count 1:2 Is Major Partner of 1",
                        "Major Partners Count 2:Mutual Major Trade Partners",
                        "Major Partners Count 2:1 Is Major Partner of 2", "Major Partners Count 2:2 Is Major Partner of 1",
                        "Trade/GDP 1:Trade/GDP 2", "Democracy 1:Democracy 2")

coefplot_f_d_l_count <- coef_fun(termname=termname_f_d_l_count, model=fitsur.f.d.l.count)

ggsave(file = "AppendixMaterials/DYLDVModelcount/coefplot_f_d_l_count.pdf", coefplot_f_d_l_count, width=10, height=6)


## Marginal Impact Plots ------------------------------------------------------

margin_plot_bothmajor_f_d_l_count <- margin_fun(data=final.data.dyad.count, model=fitsur.f.d.l.count,para="major.sym")
margin_plot_selfmajor_f_d_l_count <- margin_fun(data=final.data.dyad.count, model=fitsur.f.d.l.count,para="major.asym1")
margin_plot_targetmajor_f_d_l_count <- margin_fun(data=final.data.dyad.count, model=fitsur.f.d.l.count,para="major.asym2")
margin_plot_f_d_l_count <- grid_arrange_shared_legend(margin_plot_bothmajor_f_d_l_count+labs(title="(a) MTP both"), 
                                                    margin_plot_selfmajor_f_d_l_count +labs(title="(b) MTP self"), 
                                                    margin_plot_targetmajor_f_d_l_count+labs(title="(c) MTP target"), nrow=1)
ggsave(file = "AppendixMaterials/DYLDVModelcount/marginplot_f_d_l_count.pdf", margin_plot_f_d_l_count, width=12, height=6)

margin_plot_kof_f_d_l_count <- margin_fun2(model=fitsur.f.d.l.count, data=final.data.dyad.count)
ggsave(file = "AppendixMaterials/DYLDVModelcount/marginplot_kof_f_d_l_count.pdf", margin_plot_kof_f_d_l_count, width=12, height=6)

margin_plot_allin1_f_d_l_count <- grid.arrange(margin_plot_f_d_l_count,
                                             margin_plot_kof_f_d_l_count+labs(title="(d) KOF"),
                                             nrow=2)
ggsave(file = "AppendixMaterials/DYLDVModelcount/marginplot_allin1_f_d_l_count.pdf", margin_plot_allin1_f_d_l_count, width=12, height=10)




## Prediction Plot ------------------------------------------------------------

predict_plot_f_d_l_count <- predict_fun(data=final.data.dyad.count, model=fitsur.f.d.l.count)
ggsave(file = "AppendixMaterials/DYLDVModelcount/predictplot_f_d_l_count.pdf", predict_plot_f_d_l_count, width=10, height=6)

## Table
latex_fun(model=fitsur.f.d.l.count, modelname="DYLDVModelcount", data=final.data.dyad.count, formfake=coop.f.d.l.count, termname=termname_f_d_l_count, labelname=c("Cooperation", "Conflict"))


################################
## 8.b DY LDV KOF replaced with major partners count (both import and export)
################################
dir.create("./AppendixMaterials/DYLDVModelcountboth")

## coef plot -----------------------------------------------------------------
## term names for DY LDV model
termname_f_d_l_countboth <- c("Coop Level","Conf Level","Major Partners Count 1","Major Partners Count 2",
                          "Mutual Major Trade Partners", "1 Is Major Partner of 2", "2 Is Major Partner of 1", 
                          "Trade/GDP 1", "Trade/GDP 2", "Trade Balance",
                          "CINC Ratio", "log Distance", "log Population 1",
                          "log Population 2", "log GDP 1", "log GDP 2",
                          "Democracy 1", "Democracy 2",
                          "Alliance","Contiguity",
                          "Major Partners Count 1:Mutual Major Trade Partners",
                          "Major Partners Count 1:1 Is Major Partner of 2", "Major Partners Count 1:2 Is Major Partner of 1",
                          "Major Partners Count 2:Mutual Major Trade Partners",
                          "Major Partners Count 2:1 Is Major Partner of 2", "Major Partners Count 2:2 Is Major Partner of 1",
                          "Trade/GDP 1:Trade/GDP 2", "Democracy 1:Democracy 2")

coefplot_f_d_l_countboth <- coef_fun(termname=termname_f_d_l_countboth, model=fitsur.f.d.l.countboth)

ggsave(file = "AppendixMaterials/DYLDVModelcountboth/coefplot_f_d_l_countboth.pdf", coefplot_f_d_l_countboth, width=10, height=6)


## Marginal Impact Plots ------------------------------------------------------

margin_plot_bothmajor_f_d_l_countboth <- margin_fun(data=final.data.dyad.countboth, model=fitsur.f.d.l.countboth,para="major.sym")
margin_plot_selfmajor_f_d_l_countboth <- margin_fun(data=final.data.dyad.countboth, model=fitsur.f.d.l.countboth,para="major.asym1")
margin_plot_targetmajor_f_d_l_countboth <- margin_fun(data=final.data.dyad.countboth, model=fitsur.f.d.l.countboth,para="major.asym2")
margin_plot_f_d_l_countboth <- grid_arrange_shared_legend(margin_plot_bothmajor_f_d_l_countboth+labs(title="(a) MTP both"), 
                                                      margin_plot_selfmajor_f_d_l_countboth +labs(title="(b) MTP self"), 
                                                      margin_plot_targetmajor_f_d_l_countboth+labs(title="(c) MTP target"), nrow=1)
ggsave(file = "AppendixMaterials/DYLDVModelcountboth/marginplot_f_d_l_countboth.pdf", margin_plot_f_d_l_countboth, width=12, height=6)

margin_plot_kof_f_d_l_countboth <- margin_fun2(model=fitsur.f.d.l.countboth, data=final.data.dyad.countboth)
ggsave(file = "AppendixMaterials/DYLDVModelcountboth/marginplot_kof_f_d_l_countboth.pdf", margin_plot_kof_f_d_l_countboth, width=12, height=6)

margin_plot_allin1_f_d_l_countboth <- grid.arrange(margin_plot_f_d_l_countboth,
                                               margin_plot_kof_f_d_l_countboth+labs(title="(d) KOF"),
                                               nrow=2)
ggsave(file = "AppendixMaterials/DYLDVModelcountboth/marginplot_allin1_f_d_l_countboth.pdf", margin_plot_allin1_f_d_l_countboth, width=12, height=10)



## Prediction Plot ------------------------------------------------------------

predict_plot_f_d_l_countboth <- predict_fun(data=final.data.dyad.countboth, model=fitsur.f.d.l.countboth)
ggsave(file = "AppendixMaterials/DYLDVModelcountboth/predictplot_f_d_l_countboth.pdf", predict_plot_f_d_l_countboth, width=10, height=6)

## Table
latex_fun(model=fitsur.f.d.l.countboth, modelname="DYLDVModelcountboth", data=final.data.dyad.countboth, formfake=coop.f.d.l.countboth, termname=termname_f_d_l_countboth, labelname=c("Cooperation", "Conflict"))


################################
## 9. DY LDV with 4 lagged coop and conflict 
################################
dir.create("./AppendixMaterials/DYLDVModel4iv")

## coef plot -----------------------------------------------------------------
## term names for DY LDV model
termname_f_d_l_4iv <- c("Coop Level 1","Conf Level 1",
                        "Coop Level 2","Conf Level 2",
                        "KOF 1","KOF 2",
                        "Mutual Major Trade Partners", "1 Is Major Partner of 2", "2 Is Major Partner of 1", 
                        "Trade/GDP 1", "Trade/GDP 2", "Trade Balance",
                        "CINC Ratio", "log Distance", "log Population 1",
                        "log Population 2", "log GDP 1", "log GDP 2",
                        "Democracy 1", "Democracy 2",
                        "Alliance","Contiguity",
                        "KOF 1:Mutual Major Trade Partners",
                        "KOF 1:1 Is Major Partner of 2", "KOF 1:2 Is Major Partner of 1",
                        "KOF 2:Mutual Major Trade Partners",
                        "KOF 2:1 Is Major Partner of 2", "KOF 2:2 Is Major Partner of 1",
                        "Trade/GDP 1:Trade/GDP 2", "Democracy 1:Democracy 2")

coefplot_f_d_l_4iv <- coef_fun(termname=termname_f_d_l_4iv, model=fitsur.f.d.l.4iv)

ggsave(file = "AppendixMaterials/DYLDVModel4iv/coefplot_f_d_l_4iv.pdf", coefplot_f_d_l_4iv, width=10, height=6)


## Marginal Impact Plots ------------------------------------------------------

margin_plot_bothmajor_f_d_l_4iv <- margin_fun(data=final.data.dyad, model=fitsur.f.d.l.4iv,para="major.sym")
margin_plot_selfmajor_f_d_l_4iv <- margin_fun(data=final.data.dyad, model=fitsur.f.d.l.4iv,para="major.asym1")
margin_plot_targetmajor_f_d_l_4iv <- margin_fun(data=final.data.dyad, model=fitsur.f.d.l.4iv,para="major.asym2")
margin_plot_f_d_l_4iv <- grid_arrange_shared_legend(margin_plot_bothmajor_f_d_l_4iv+labs(title="(a) MTP both"), 
                                                      margin_plot_selfmajor_f_d_l_4iv +labs(title="(b) MTP self"), 
                                                      margin_plot_targetmajor_f_d_l_4iv+labs(title="(c) MTP target"), nrow=1)
ggsave(file = "AppendixMaterials/DYLDVModel4iv/marginplot_f_d_l_4iv.pdf", margin_plot_f_d_l_4iv, width=12, height=6)


margin_plot_kof_f_d_l_4iv <- margin_fun2(model=fitsur.f.d.l.4iv, data=final.data.dyad)
ggsave(file = "AppendixMaterials/DYLDVModel4iv/marginplot_kof_f_d_l_4iv.pdf", margin_plot_kof_f_d_l_4iv, width=12, height=6)

margin_plot_allin1_f_d_l_4iv <- grid.arrange(margin_plot_f_d_l_4iv,
                                                   margin_plot_kof_f_d_l_4iv+labs(title="(d) KOF"),
                                                   nrow=2)
ggsave(file = "AppendixMaterials/DYLDVModel4iv/marginplot_allin1_f_d_l_4iv.pdf", margin_plot_allin1_f_d_l_4iv, width=12, height=10)




## Prediction Plot ------------------------------------------------------------

predict_plot_f_d_l_4iv <- predict_fun(data=final.data.dyad, model=fitsur.f.d.l.4iv)
ggsave(file = "AppendixMaterials/DYLDVModel4iv/predictplot_f_d_l_4iv.pdf", predict_plot_f_d_l_4iv, width=10, height=6)

## Table
latex_fun(model=fitsur.f.d.l.4iv, modelname="DYLDVModel4iv", data=final.data.dyad, formfake=coop.f.d.l.4iv, termname=termname_f_d_l_4iv, labelname=c("Cooperation", "Conflict"))

################################
## 10. DY LDV KOF replaced with dep
################################
dir.create("./AppendixMaterials/DYLDVModeldep")

## coef plot -----------------------------------------------------------------
## term names for DY LDV model
termname_f_d_l_dep <- c("Coop Level","Conf Level","Trade/GDP 1","Trade/GDP 2",
                          "Mutual Major Trade Partners", "1 Is Major Partner of 2", "2 Is Major Partner of 1", 
                          "Trade Balance",
                          "CINC Ratio", "log Distance", "log Population 1",
                          "log Population 2", "log GDP 1", "log GDP 2",
                          "Democracy 1", "Democracy 2",
                          "Alliance","Contiguity",
                          "Trade/GDP 1:Mutual Major Trade Partners",
                          "Trade/GDP 1:1 Is Major Partner of 2", "Trade/GDP 1:2 Is Major Partner of 1",
                          "Trade/GDP 2:Mutual Major Trade Partners",
                          "Trade/GDP 2:1 Is Major Partner of 2", "Trade/GDP 2:2 Is Major Partner of 1",
                           "Democracy 1:Democracy 2")

coefplot_f_d_l_dep <- coef_fun(termname=termname_f_d_l_dep, model=fitsur.f.d.l.dep)

ggsave(file = "AppendixMaterials/DYLDVModeldep/coefplot_f_d_l_dep.pdf", coefplot_f_d_l_dep, width=10, height=6)


## Marginal Impact Plots ------------------------------------------------------

margin_plot_bothmajor_f_d_l_dep <- margin_fun(data=final.data.dyad.dep, model=fitsur.f.d.l.dep,para="major.sym")
margin_plot_selfmajor_f_d_l_dep <- margin_fun(data=final.data.dyad.dep, model=fitsur.f.d.l.dep,para="major.asym1")
margin_plot_targetmajor_f_d_l_dep <- margin_fun(data=final.data.dyad.dep, model=fitsur.f.d.l.dep,para="major.asym2")
margin_plot_f_d_l_dep <- grid_arrange_shared_legend(margin_plot_bothmajor_f_d_l_dep+labs(title="(a) MTP both"), 
                                                      margin_plot_selfmajor_f_d_l_dep +labs(title="(b) MTP self"), 
                                                      margin_plot_targetmajor_f_d_l_dep+labs(title="(c) MTP target"), nrow=1)
ggsave(file = "AppendixMaterials/DYLDVModeldep/marginplot_f_d_l_dep.pdf", margin_plot_f_d_l_dep, width=12, height=6)


margin_plot_kof_f_d_l_dep <- margin_fun2(model=fitsur.f.d.l.dep, data=final.data.dyad.dep)
ggsave(file = "AppendixMaterials/DYLDVModeldep/marginplot_kof_f_d_l_dep.pdf", margin_plot_kof_f_d_l_dep, width=12, height=6)

margin_plot_allin1_f_d_l_dep <- grid.arrange(margin_plot_f_d_l_dep,
                                             margin_plot_kof_f_d_l_dep+labs(title="(d) KOF"),
                                             nrow=2)
ggsave(file = "AppendixMaterials/DYLDVModeldep/marginplot_allin1_f_d_l_dep.pdf", margin_plot_allin1_f_d_l_dep, width=12, height=10)



## Prediction Plot ------------------------------------------------------------

predict_plot_f_d_l_dep <- predict_fun(data=final.data.dyad.dep, model=fitsur.f.d.l.dep)
ggsave(file = "AppendixMaterials/DYLDVModeldep/predictplot_f_d_l_dep.pdf", predict_plot_f_d_l_dep, width=10, height=6)

## Table
latex_fun(model=fitsur.f.d.l.dep, modelname="DYLDVModeldep", data=final.data.dyad.dep, formfake=coop.f.d.l.dep, termname=termname_f_d_l_dep, labelname=c("Cooperation", "Conflict"))

################################
## 11. DY LDV with third party controls
################################
dir.create("./AppendixMaterials/DYLDVModel3rdevents")

## coef plot -----------------------------------------------------------------
## term names for DY LDV model
termname_f_d_l_3rdevents <- c("Coop Level","Conf Level",
                        "Coop with Others 1","Conf with Others 1",
                        "Coop with Others 2","Conf with Others 2",
                        "KOF 1","KOF 2",
                        "Mutual Major Trade Partners", "1 Is Major Partner of 2", "2 Is Major Partner of 1", 
                        "Trade/GDP 1", "Trade/GDP 2", "Trade Balance",
                        "CINC Ratio", "log Distance", "log Population 1",
                        "log Population 2", "log GDP 1", "log GDP 2",
                        "Democracy 1", "Democracy 2",
                        "Alliance","Contiguity",
                        "KOF 1:Mutual Major Trade Partners",
                        "KOF 1:1 Is Major Partner of 2", "KOF 1:2 Is Major Partner of 1",
                        "KOF 2:Mutual Major Trade Partners",
                        "KOF 2:1 Is Major Partner of 2", "KOF 2:2 Is Major Partner of 1",
                        "Trade/GDP 1:Trade/GDP 2", "Democracy 1:Democracy 2")

coefplot_f_d_l_3rdevents <- coef_fun(termname=termname_f_d_l_3rdevents, model=fitsur.f.d.l.3rdevents)

ggsave(file = "AppendixMaterials/DYLDVModel3rdevents/coefplot_f_d_l_3rdevents.pdf", coefplot_f_d_l_3rdevents, width=10, height=6)


## Marginal Impact Plots ------------------------------------------------------

margin_plot_bothmajor_f_d_l_3rdevents <- margin_fun(data=final.data.dyad, model=fitsur.f.d.l.3rdevents,para="major.sym")
margin_plot_selfmajor_f_d_l_3rdevents <- margin_fun(data=final.data.dyad, model=fitsur.f.d.l.3rdevents,para="major.asym1")
margin_plot_targetmajor_f_d_l_3rdevents <- margin_fun(data=final.data.dyad, model=fitsur.f.d.l.3rdevents,para="major.asym2")
margin_plot_f_d_l_3rdevents <- grid_arrange_shared_legend(margin_plot_bothmajor_f_d_l_3rdevents+labs(title="(a) MTP both"), 
                                                    margin_plot_selfmajor_f_d_l_3rdevents +labs(title="(b) MTP self"), 
                                                    margin_plot_targetmajor_f_d_l_3rdevents+labs(title="(c) MTP target"), nrow=1)
ggsave(file = "AppendixMaterials/DYLDVModel3rdevents/marginplot_f_d_l_3rdevents.pdf", margin_plot_f_d_l_3rdevents, width=12, height=6)

margin_plot_kof_f_d_l_3rdevents <- margin_fun2(model=fitsur.f.d.l.3rdevents, data=final.data.dyad)
ggsave(file = "AppendixMaterials/DYLDVModel3rdevents/marginplot_kof_f_d_l_3rdevents.pdf", margin_plot_kof_f_d_l_3rdevents, width=12, height=6)

margin_plot_allin1_f_d_l_3rdevents <- grid.arrange(margin_plot_f_d_l_3rdevents,
                                             margin_plot_kof_f_d_l_3rdevents+labs(title="(d) KOF"),
                                             nrow=2)
ggsave(file = "AppendixMaterials/DYLDVModel3rdevents/marginplot_allin1_f_d_l_3rdevents.pdf", margin_plot_allin1_f_d_l_3rdevents, width=12, height=10)



## Prediction Plot ------------------------------------------------------------

predict_plot_f_d_l_3rdevents <- predict_fun(data=final.data.dyad, model=fitsur.f.d.l.3rdevents)
ggsave(file = "AppendixMaterials/DYLDVModel3rdevents/predictplot_f_d_l_3rdevents.pdf", predict_plot_f_d_l_3rdevents, width=10, height=6)

## Table
latex_fun(model=fitsur.f.d.l.3rdevents, modelname="DYLDVModel3rdevents", data=final.data.dyad, formfake=coop.f.d.l.3rdevents, termname=termname_f_d_l_3rdevents, labelname=c("Cooperation", "Conflict"))


#######################
## 12. Exclude zero interactions

dir.create("./AppendixMaterials/DYLDVwozeroModel")

## coef plot -----------------------------------------------------------------
## term names for DY LDV model
termname_f_d_l_wozero <- c("Coop Level","Conf Level","KOF 1","KOF 2",
                    "Mutual Major Trade Partners", "1 Is Major Partner of 2", "2 Is Major Partner of 1", 
                    "Trade/GDP 1", "Trade/GDP 2", "Trade Balance",
                    "CINC Ratio", "log Distance", "log Population 1",
                    "log Population 2", "log GDP 1", "log GDP 2",
                    "Democracy 1", "Democracy 2",
                    "Alliance","Contiguity",
                    "KOF 1:Mutual Major Trade Partners",
                    "KOF 1:1 Is Major Partner of 2", "KOF 1:2 Is Major Partner of 1",
                    "KOF 2:Mutual Major Trade Partners",
                    "KOF 2:1 Is Major Partner of 2", "KOF 2:2 Is Major Partner of 1",
                    "Trade/GDP 1:Trade/GDP 2", "Democracy 1:Democracy 2")

coefplot_f_d_l_wozero <- coef_fun(termname=termname_f_d_l_wozero, model=fitsur.f.d.l.wozero)

ggsave(file = "AppendixMaterials/DYLDVwozeroModel/coefplot_f_d_l_wozero.pdf", coefplot_f_d_l_wozero, width=10, height=6)


## Marginal Impact Plots ------------------------------------------------------

margin_plot_bothmajor_f_d_l_wozero <- margin_fun(model=fitsur.f.d.l.wozero,para="major.sym",data=final.data.dyad %>% filter(coopTotal!=0|hostTotal!=0))
margin_plot_selfmajor_f_d_l_wozero <- margin_fun(model=fitsur.f.d.l.wozero,para="major.asym1",data=final.data.dyad %>% filter(coopTotal!=0|hostTotal!=0))
margin_plot_targetmajor_f_d_l_wozero <- margin_fun(model=fitsur.f.d.l.wozero,para="major.asym2",data=final.data.dyad %>% filter(coopTotal!=0|hostTotal!=0))
margin_plot_f_d_l_wozero <- grid_arrange_shared_legend(margin_plot_bothmajor_f_d_l_wozero+labs(title="(a) MTP both"), 
                                                margin_plot_selfmajor_f_d_l_wozero +labs(title="(b) MTP self"), 
                                                margin_plot_targetmajor_f_d_l_wozero+labs(title="(c) MTP target"), nrow=1)
ggsave(file = "AppendixMaterials/DYLDVwozeroModel/marginplot_f_d_l_wozero.pdf", margin_plot_f_d_l_wozero, width=12, height=6)

margin_plot_kof_f_d_l_wozero <- margin_fun2(model=fitsur.f.d.l.wozero,data=final.data.dyad %>% filter(coopTotal!=0|hostTotal!=0))
ggsave(file = "AppendixMaterials/DYLDVwozeroModel/marginplot_kof_f_d_l_wozero.pdf", margin_plot_kof_f_d_l_wozero, width=12, height=6)

margin_plot_allin1_f_d_l_wozero <- grid.arrange(margin_plot_f_d_l_wozero,
                                         margin_plot_kof_f_d_l_wozero+labs(title="(d) KOF"),
                                         nrow=2)
ggsave(file = "AppendixMaterials/DYLDVwozeroModel/marginplot_allin1_f_d_l_wozero.pdf", margin_plot_allin1_f_d_l_wozero, width=12, height=10)


## Prediction Plot ------------------------------------------------------------

predict_plot_f_d_l_wozero <- predict_fun(model=fitsur.f.d.l.wozero,data=final.data.dyad %>% filter(coopTotal!=0|hostTotal!=0))
ggsave(file = "AppendixMaterials/DYLDVwozeroModel/predictplot_f_d_l_wozero.pdf", predict_plot_f_d_l_wozero, width=10, height=6)


## Main Table: DY LDV ---------------------------------------------------------
latex_fun(model=fitsur.f.d.l.wozero, modelname="DYLDVwozeroModel", modelcap="lagged events model", 
          data=final.data.dyad %>% filter(coopTotal!=0|hostTotal!=0), 
          formfake=coop.f.d.l, termname=termname_f_d_l_wozero, 
          labelname=c("Cooperation", "Conflict"))


###########################
## 13. DY KOF FE
###########################
dir.create("./AppendixMaterials/DYFEKOFModel")

## Coef Plot
termname_f_d_kof <- c("KOF 1","KOF 2","Mutual Major Trade Partners",
                  "1 Is Major Partner of 2", "2 Is Major Partner of 1", "Trade/GDP 1", "Trade/GDP 2",
                  "Trade Balance",
                  "CINC Ratio", "log Distance", "log Population 1",
                  "log Population 2", "log GDP PC 1", "log GDP PC 2",
                  "Democracy 1", "Democracy 2", "KOF yearly average","Alliance","Contiguity", "KOF 1:Mutual Major Trade Partners",
                  "KOF 1:1 Is Major Partner of 2", "KOF 1:2 Is Major Partner of 1",
                  "KOF 2:Mutual Major Trade Partners",
                  "KOF 2:1 Is Major Partner of 2", "KOF 2:2 Is Major Partner of 1",
                  "Trade/GDP 1:Trade/GDP 2", "Democracy 1:Democracy 2")

coefplot_f_d_kof <- coef_fun(termname=termname_f_d_kof, model=fitsur.f.d.kof)

ggsave(file = "AppendixMaterials/DYFEKOFModel/coefplot_f_d_kof.pdf", coefplot_f_d_kof, width=10, height=6)


## Marginal Effects plot
margin_plot_bothmajor_f_d_kof <- margin_fun(model=fitsur.f.d.kof,para="major.sym")
margin_plot_selfmajor_f_d_kof <- margin_fun(model=fitsur.f.d.kof,para="major.asym1")
margin_plot_targetmajor_f_d_kof <- margin_fun(model=fitsur.f.d.kof,para="major.asym2")
margin_plot_f_d_kof <- grid_arrange_shared_legend(margin_plot_bothmajor_f_d_kof+labs(title="(a) MTP both"), 
                                              margin_plot_selfmajor_f_d_kof+labs(title="(b) MTP self"), 
                                              margin_plot_targetmajor_f_d_kof+labs(title="(c) MTP target"), nrow=1)
ggsave(file = "AppendixMaterials/DYFEKOFModel/marginplot_f_d_kof.pdf", margin_plot_f_d_kof, width=12, height=6)

margin_plot_kof_f_d_kof <- margin_fun2(model=fitsur.f.d.kof)
ggsave(file = "AppendixMaterials/DYFEKOFModel/marginplot_kof_f_d_kof.pdf", margin_plot_kof_f_d_kof, width=12, height=6)

margin_plot_allin1_f_d_kof <- grid.arrange(margin_plot_f_d_kof,
                                       margin_plot_kof_f_d_kof+labs(title="(d) KOF"),
                                       nrow=2)
ggsave(file = "AppendixMaterials/DYFEKOFModel/marginplot_allin1_f_d_kof.pdf", margin_plot_allin1_f_d_kof, width=12, height=10)


## Predict Plot
predict_plot_f_d_kof <- predict_fun(model=fitsur.f.d.kof, 
                                    data=final.data.dyad)
ggsave(file = "AppendixMaterials/DYFEKOFModel/predictplot_f_d_kof.pdf", predict_plot_f_d_kof, width=10, height=6)


## Table
latex_fun(model=fitsur.f.d.kof, modelname="DYFEKOFModel", 
          modelcap="model with KOF yearly average", 
          data=final.data.dyad, formfake=coop.f.d.kof, 
          termname=termname_f_d_kof, labelname=c("Cooperation", "Conflict"))

